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Abstract 

Combining recent techniques giving non-perturbative re-summed estimates 
of the damping and conservative parts of the two-body dynamics, we describe 
the transition between the adiabatic phase and the plunge, in coalescing bi- 
nary black holes with comparable masses moving on quasi-circular orbits. We 
give initial dynamical data for numerical relativity investigations, with a frac- 
tion of an orbit left, and provide, for data analysis purposes, an estimate 
of the gravitational wave-form emitted throughout the inspiral, plunge and 
coalescence phases. 
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I. INTRODUCTION 



The most promising candidate sources for ground based interferometric gravitational- 
wave (GW) detectors such as LIGO and VIRGO are binary systems made of massive (stellar) 
black holes Such binary black holes (with individual masses in the range, say, 3— I5M0) 

pose special problems @, f^. 

Let us recall that gravitational radiation damping is efficient at circularizing such bi- 
nary systems, and then drives, for a long time, a slow inspiraling quasi- circular motion of 
the binary system. This quasi-circular "adiabatic inspiral phase" is expected to terminate 
abruptly, and to change to some type of "plunge phase" (leading to final coalescence) when 
the binary orbit shrinks down to the Last Stable (circular) Orbit (LSO) defined by the conser- 
vative part of the nonlinear relativistic force law between two bodies. [In the test-mass limit, 
the full nonlinear relativistic force law corresponds to geodesic motion in a Schwarzschild 
spacetime, and exhibits, as is well known, an LSO located at i? = 6GM. One expects that 
a comparable-mass system will still exhibit such an LSO; see below.] Now, the signal to 
noise ratio (in an initial LIGO detector) for inspiral signals from comparable-mass black 
hole binaries reaches a maximum for M ~ 28Mq, which corresponds to a GW frequency for 
the waves emitted at the Last Stable (circular) Orbit (LSO) equal to /gw* — 170IIz, a value 
which is (not accidentally) very close to the location /^ct — 167IIz (for initial LIGO) of the 
minimum of the characteristic detector noise amplitude /i„(/) = a/ / Sn{f) (see Fig. |^ of 
0). Therefore the first detections will probably concern massive systems with M ~ SOMq. 
Moreover, Ref. |^ has shown that when the total mass M = mi + m2 lies in the range 
5 — 40Mq the proximity (within a factor ~ 2) between the observationally most important 
frequencies^] fdet and the GW frequency at the LSO, > calling both for an especially 
careful treatment of the Fourier transform of the emitted waveform, and for an improved 
knowledge of the transition between the inspiral phase and the plunge phase. 

The present paper will attempt to improve our knowledge of the transition between inspi- 
ral and plunge by combining two, recently proposed, non-perturbative techniques: Refs. 



We neglect here the very small difference between the optimal frequency fdet for generic broad- 
band bursts, and the optimal frequency fp for inspiral signals (see |Q). 
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and 0. Let us first recall that, a few years ago, Will and collaborators [Q, tried to 
attack the problem of the late-time evolution of compact binaries (including the transition 
from inspiral to plunge) by a direct use of the D amour- Deruelle |TO|-p!^ equations of motion. 
These equations of motion are given in the form of a perturbative expansion in powers of 
a small parameter e = v/c ("post-Newtonian", or, in short, PN, expansion). In Ref. a 
direct integration of these perturbative equations of motion (using the method of osculating 
elements) was used, while, in Ref. |^ it was proposed to improve the straightforward pertur- 
bative approach by using "hybrid" equations of motion. The "hybrid" approach is a partial 
re-summation approach in which the perturbative terms in the equations of motion which 
survive in the test mass limit {u = mi 7712/ (mi +7712)^ 0) are replaced by the known, exact 
"Schwarzschild terms", while the zz-dependent terms are left as a perturbative expansion. 
Both the robustness [|13|, |T^ and the consistency of the hybrid approach of have 
been questioned. [In particular, it was pointed out in Ref. that, in this approach, the 
supposedly small "//-corrections" represent, in several cases, a very large (larger than 100%) 
modification of the corresponding //-independent terms.] Another sign of the unreliability of 
the hybrid approach is the fact that the recent study pr5|Jl6|] of the location of the LSO at 
the third post-Newtonian (3PN) accuracy has qualitatively confirmed the 2PN-level results 
of the non-perturbative techniques to be discussed below (namely that the LSO is "lower 
than 6GM"), thereby casting doubt on the most striking prediction of the hybrid approach 
(an LSO "higher than 6GM", i.e. with a lower orbital frequency). 

By contrast with the perturbative approach of and the partially re-summed approach 
of 0, the present paper will rely on the systematic use of non-perturbative re-summation 
techniques. The basic philosophy underlying our approach is the following. We are interested 
in understanding, in quantitative detail, the combined infiuence on the inspiral plunge 
transition of radiation reaction and of non-linear effects in the force law for comparable-mass 
binary systems (i.e. for systems for which u = nii m^j irrix + m^)^ is around^ 1/4). 

At present there exists no method for deriving, from first principles, non-perturbative 



^ Note that because z/, considered as a function of the ratio mi/m2, reaches its maximum fmax = 
1/4 for mi/m2 = 1 it stays numerically near 1/4 even for mass ratios quite different from 1. E.g., 
even for mxjm^ = 3, 4i/ = 0.75. 
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expressions for the two-body equations of motion, especially in the case of interest where 

is not a small parameter. As a substitute we shall combine two different re-summation 
techniques that have been recently introduced to deal with two separate aspects of the 
problem we wish to tackle. 

The first re-summation technique, introduced in 0, allows one to get a non-perturbative, 
z/-dependent, estimate of the rate of loss of angular momentum (under gravitational damping) 
in quasi-circular, comparable-mass binaries. The idea of |0 is three-pronged: (i) to work 
with an invariant function of an invariant argument, F{v), (ii) to inject some plausible 
information about the meromorphic structure of this function, and, finally, (iii) to use Fade 
approximants to estimate F{v) from the first few known terms in the perturbative (PN) 
expansion of F{v). The second re-summation technique, introduced in 0, allows one to 
derive a non-perturbative, i/-dependent, estimate of the (conservative part of the) nonlinear 
force law determining the motion of comparable binaries. The idea of [0] is to map the real 
two-body problem on a simpler effective one-body problem, i.e. on the problem of the motion 
of a particle of mass fi = 11111712/ {mi + in some "effective" background metric gj^ix^). 
The possibility (and uniqueness, given some natural requirement) of such a mapping, real 

effective, was proven at the 2PN level in 0. The extension of this approach at the 3PN 
level has been recently discussed |1I5|JT6[| on the basis of the 3PN dynamics recently derived 
by Jaranowski and Schafer [0. At the 2PN level the i/-dependent terms in the effective 
metric were found to be numerically so small (around the LSO) that the need for a further 
(Pade-type) re-summing of the effective metric coefficients did not arise. [However, note that 
Ref. [0 has introduced, at 2PN and 3PN, the further idea of a specific, Pade-improvement 
of gj^ix^).] In this paper we shall show how one can combine the methods of and [0 
to derive a full force law (including radiation reaction) describing the quasi-circular motion 
of comparable-mass binaries. Our approach is intended to apply to any value of u, but is 
restricted to considering quasi- circular motions, where the radial velocity R is much smaller 
than the circular one R(f. As we shall see, we shall consistently check that the condition 
R <^ Rip holds true not only during the adiabatic inspiral, but also during the transition to 
the plunge, and even during most of the plunge. 

We apply our method, in this paper, to deriving two sorts of results which are of di- 
rect interest to the ongoing effort to detect gravitational waves. First, we shall give initial 
dynamical data (i.e. initial positions and momenta) for binary black holes that have just 
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started their plunge motion. The idea here is that numerical relativity will probably not 
be able, before quite a few years, to accurately evolve binary systems over many (or even 
~ 10) orbits. This is why we propose a method for computing accurate initial dynamical 
data at a moment so late in the evolution that there remains (when 4i/ ~ 1) less than one 
orbit to evolve. [In the equal-mass case, z/ = 1/4, we shall compute data ~ 0.6 orbit before 
"coalescence" .] Our contention (whose robustness we shall try to establish) is that suitably 
re-summed versions of analytic (PN) results allow one to push the evolution that far. [We 
shall use here 2.5PN-accurate information for angular momentum loss and 2PN- accurate 
information for the conservative force law. However, as shown in [§] and |1T5| , p!6| our method 
can be pushed to higher accuracy when the correspondingly needed PN results become un- 
ambiguously known.] Note that this attitude is opposite to the one taken in [^ in which 
it was assumed that "there is little hope, via PN Pade approximants, to evolve" a binary 
system up to the moment where it can provide initial data for the final coalescence. Let us, 
however, immediately add that the present paper is still incomplete, in that we give only 
dynamical data (q^i, Pi5 P2) but we do not solve the remaining problem of constructing 
the initial gravitational data {gij{x), Kij{x)) determined (in principle) by {qa,Pa) (given, 
say, some no-incoming-radiation condition). We shall leave this (important) issue to future 
work. 

The second aim of this work is to provide, for data analysis purposes, some estimate of 
the complete waveform emitted by the coalescence of two black holes (with negligible spins). 
We do not claim that this part of the work will be as accurate as the first one. The idea 
here is to provide a (hopefully ~ 10% accurate) guess of the complete waveform, with its 
transition from an inspiral phase to a plunge one, followed by a coalescence ending in a 
stationary final state. In view of the recent realization [^ of the crucial importance of the 
details of the transition to the plunge for the construction of faithful GW templates (for 
massive binaries with 5Mq <M< 4OM0) even an approximate knowledge of the complete 
waveform will be a valuable information for data analysis (e.g. to test the accuracy of present 
templates, and/or to propose more accurate or, at least, more robust, templates). 

While preparing this work for publication, we learned of the existence of an independent 
work of Ori and Thorne |jl8[ which deals with the transition between the inspiral and the 
plunge in the test mass limit ( — > 0). 
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II. CONSERVATIVE PART OF THE TWO-BODY FORCE LAW 



In this section, we recall the non-perturbative construction of the (conservative) two- 
body force law given in Ref. 0. There it was shown that the conservative part (i.e. 
without radiation damping) of the dynamics of a binary system, represented in ADM 
phase-space coordinates (qf°^, g^°'^,pf°^,p^°'^), could be mapped (at the 2PN level), 
via the combination of an energy map, Ses = /(^reai), and a canonical transformation, 
^gADM^pADM^ _^ (^aiPa)? = 1; 2, iuto the simpler dynamics of the geodesic motion of a 
particle of mass fi = mi 1712/ {mi + 7712) in some effective background geometry gj^lx): 

dsls = gf^x^) dx^ dx" = -A{R) + B{R) dR^ + C{R) R^ {dO^ + sin^ 6 d^^) . (2.1) 

[See [0 for the generalization of this approach to the 3PN level.] Here the coordinates 
{R, 6, ip) are polar coordinates in the effective problem (describing the relative motion). They 
are related in the standard way {Q^ = R sin 6 cosip, = R sin 6 sinp, = R cos 6) to 
the (relative) effective Cartesian coordinates Q = Qi — 0,21 where qi and ^2 are the effective 
coordinates of each body. One works in the center-of-mass frame of the binary system, i.e. 
Pi ^ P2 — ^ — Pi^^ + P2^^. The canonical conjugate of the relative position Q is the 
relative momentum P = = —p2- In most of this paper we shall work with the effective 
phase-space coordinates {Q, P) [or rather with scaled versions of their polar f\ counterparts 
{R, 9, ip; Pr, Pg, -P^)]. We shall only discuss at the end how to construct the more physically 
relevant ADM phase space coordinates {Qa^^yPa^^) from {Q,P). 

In absence of damping (to be added later), the evolution (with respect to the real ADM 
time coordinate treai) of (Q, P) is given by Hamilton's equations 

dQ^ 9/73°"" (Q.^) 



dtj-eal dPi 



0, (2.2) 



where the real (i.e. giving the treai-evolution, and the real two-body energy) improved (i.e. 
representing a non-perturbative re-summed estimate of the real PN Hamiltonian) Hamilto- 
nian reads 



^ Note that we have the usual relations, such as, Pr = with = Q^/R, and Pip = Py 

QyPx. 
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improved 
real 



(2.4) 



and 



H,siQ,P)^f^c\ A{Q) 



1 + 



(n-P) 



in X py 



(2.5) 



Here Q = •>/ 5^ Qi — R, — Q'^/Q is the unit vector in the radial direction, and the 
scalar and vector products are performed as in Euclidean space. Henceforth, we shall pose 
t = treai, H = and use the following notation: 



M = mi + m2 , 



mi m2 
M 



_ /J' _ rrii 7712 
M (mi + m2)^ 



(2.6) 



In polar coordinates, restricting ourselves to planar motion in the equatorial plane 9 
7r/2 and to the Schwarzschild gauge {C{Q) — 1), we get the equations of motion 

dR OH 



{R, Pr, p<p) , 
{R, Pr, Pif) , 



dt OPr 
dip _ dH 
'dt ~ dP'^ 



dP^ 
dt 



0, 



(2.7) 
(2.8) 

(2.9) 
2.10) 



with 



H{R,Pr,P^)^Mc^ 



\ 



l + 2u 



lA{R) 1 



p2 



^^c^B{R) fM^c^R^ 



(2.11) 



Like in any (non-degenerate) Hamiltonian system, this conservative dynamics is equivalent 
to a Lagrangian dynamics 



Li-r-^(Q, Q) ^P,Q^- H^r^^'iQ, P) , 



(2.12) 



with PiiQ) obtained by solving = dH/dPi. The Lagrangian equations of motion read: 

(2.13) 



^ Q^improved ^^improved 



dt dQi 



To ease the notation we denote L = L^^^^°^'^'^- 
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Finally, the 2PN- accurate metric coefficients A{R), B(R), Eq. (|2.1| ), (in the Schwarzschild 
gauge where C{R) = 1) read 



2GM ^ fGM\ 
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B{R)=D{R)/A{R), (2.15) 



with 

m) = ^-Gu[^) ■ (2.16) 

Note that it was recently suggested |T6[ (because of the slow convergence of the 3PN con- 
tributions) to replace the straightforward expression ( 2.14| ) by a suitably Padeed version, 



namely (at 2PN): Ap2{R) = 1 — 2u(l + uu^)^^ , where u = GM/c^ R. However, we have 
checked that this refinement has only a very minor effect on the results to be discussed below. 
The re-summed (conservative) dynamics defined by the Hamiltonian (|2.11|) contains 



a Last Stable (circular) Orbit (LSO) which is a //-deformed version of the well known 
Schwarzschild LSO. Let us recall that the radius of the LSO is obtained by imposing the 
existence of an inflection point in the effective potential H{R, Pr = 0, J') for the radial 
motion, 

— {R,P^ = 0,J) = = -^ {R, Pr = 0,J), (2.17) 



where the total angular momentum J' = P^^ stays fixed. Eq. ( p.l7| ) has a solution in R 



(for each value of u) only for some specific value of = J'^^'^'^u). In terms of the rescaled 
variables r = c^R/GM, j = cJ'/{fiGM), Q = GMi^/c^, the LSO quantities defined, in the 
equal-mass case z/ = 1/4, by the Hamiltonian ( ^.11[ ), take the following values |^ 

rLso (1/4) = 5.718 , jlso (1/4) = 3.404 , 

(^Lso (1/4) = 0.07340 , ^^j^ = -0.01501 . (2.18) 

Note that the comparable-mass LSO is slightly more inwards (both in terms of the coordinate 
R and in the sense of having a higher orbital frequency) than its corresponding rescaled test- 
mass limit: rLso(O) = 6, Jlso(O) = = 3.4641, wlso(O) = Q'^^^ = 0.068041. 
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As we shall need in the following to refer to the numerical value of u)lso{^) for arbitrary 
values of u, we have fitted the result obtained by the (rather intricate) method of Ref. [0] to 
a simple polynomial in u. We find 

Qi^soi^) ^uJo[l + uji (4z/) + UJ2 (4z/)2] , (2.19) 
ujQ = 0.0680414 , iui = 0.0693305 , cja = 0.00935142 . (2.20) 



III. INCORPORATING RADIATION REACTION EFFECTS 

We wish to augment the conservative dynamics described in the previous section by 
adding, as accurately as possible, radiation reaction effects. If we were doing it in the 
Lagrangian formalism we would write (in any coordinate system) 

This would define the additional damping force J-\'^^{Q,Q) needed in the Lagrangian for- 
malism. In particular, in polar coordinates we would write (for planar motion 9 = 7r/2): 

^|=^-W,.A.,). (3.3) 

We want to work in the Hamiltonian framework, hence coming back to the coordinates R, Pr, 
and and imposing the constraint that the usual definition P, = OL/dQ^ holds without 
corrections (which implies that the other usual relations = dH/dPi, dH/dQ^ = —dL/dQ'^ 
and Eq. ( p.l2| ) hold too) we get 

f _|£,«,P„P^,.0, (3.4, 

^_|£(^.P„P,,.0, (3.5, 

^ + |^(^, Pr, P^) = Pr, P^) , (3.6) 

^ = ^f-(P,¥.,P^,P^), (3.7) 

where the Hamiltonian damping force J^^^^{Q^ , Pj) is numerically equal to the Lagrangian 
one: ^/^"""(Q^ P,) = ^^^(Q^ = dH/dPj) . 
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A. What do we know about the radiation reaction force? 



The radiation reaction force T was computed explicitly, at lowest (Newtonian) fractional 
order, in harmonic Cartesian-like coordinates, as part of the complete 2.5PN equations of 
motion, by Damour and Deruelle [|ToH l^. An equivalent result was also derived within the 
ADM canonical formalism by Schafer pU|-pT[]. At higher post- Newtonian orders one has only 
an incomplete knowledge of the equations of motion, and one has to rely on the (assumed) 
balance between energy and angular momentum losses in the system and at infinity [p^ , ^ . 
To get an idea of the generic structure of the radiation damping (in various coordinate 
systems, and at various PN approximations) let us consider the general radiation reaction 
force written (at IPN fractional accuracy; and setting G = 1) by Iyer and Will p2 |. 

^Lag 



/i 



a(/2, v) 
P{R,v) 



a{R,v)Rn' + (3{R,v)v' 
8 M M , , 

8 M M , ^ 

Z/ — — + -D7/2 + ■ 



(3.8) 
(3.9) 
(3.10) 



5 ' i?2 ^ 

where R is the relative radius and v is the velocity. Then, using post-Newtonian expressions 
for the energy and the angular momentum flux at infinity, and assuming energy and angular 
momentum balance, they obtained at lowest (Newtonian) fractional order 



A 



5/2 = 3{l + (3) + 1 {23 + 6a -9 (3)^- 5 fSR". 

R 



^5/2 = (2 + a) v^ + {2 



a) 



M 

'r 



3(l + a)R^ 



(3.11) 
(3.12) 



See Ref. |2^ for the expressions of the IPN-accurate radiation damping terms Aj/2 and 
-B7/2 in the equations of motion (equivalent, after some reshuffling, with the Lagrangian 
contributions A'^^^, B'^j^ in Eqs. ( F^ , ( 3.10|) ). 

The coefficients a and /? that appear in Eqs. ( |3.11| ) and ( p. 12] ) are two arbitrary gauge 
parameters that cannot be fixed by the energy balance method. Iyer and Will p2[ showed 
that this gauge freedom is equivalent to shifting the (conservative) coordinate system by small 
radiative corrections. Let us notice that the gauge dependence is reduced when considering 
quasi-circular orbits. Indeed, in that case R? ~ 0, M/R ^ t>^ and Eqs. (|3.9|) , ( p. 101) become 
(considering only the 5/2PN terms which are sufficient for the point we wish to make) 



M f mV /32 „ \ ^ 32 M /M^ ^ 



/3c 



R 



(3.13) 
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Hence, in the quasi-circular case the only gauge-dependence left is in the coefficient a{R, v) 
multiplying the radial component of the damping force (oc n*). We can use this gauge 
arbitrariness to set the ratio 

1 Z' 16 \ /o , .X 

+ a , (3.14) 



(3 J . 2 V 3 

/ circ. \ 

to any value we like. For example, by choosing a = —16/3 we can set acirc. = or by 
choosing a = —10/3 we can set (a -|- /3)circ. = 0. 

Having understood the gauge dependence of the coefficient a in Eq. ( p.8| ) let us come 
back to the general structure (|3.8| ) (considered at any PN accuracy, with some (unknown) 
coefficients a and (3). The polar-coordinate version (for planar motion 9 = 7r/2) of the 
Cartesian-like Lagrangian damping force (|3.^ ) reads (|3^) , (|3.3|) with 



dR 

^Lag ^ _^Lag ^ = J-^ag _ Qy ^Lag_ (3^ ^g) 



This yields 



^Lag_^^^2^^ ^]j"s = /x (a + /?)/?. (3.17) 



The important information for our present purpose is the difference between the 
(y9-component of the damping force, which contains only (3 and is, therefore, gauge- 
independent 0, and the i?-component which contains the gauge- dependent combination a+f3. 
Let us note, in particular, the expression of the ratio 

^J;"^ fa \ R , , 

- ' ^5 + 1 (3.18) 



In the following we shall be interested in quasi-circular motions with R <^ R(f. [We shall 
see that this condition remains satisfied even during part of the plunge phase.] As we see 
from Eq. (|3.18|), for such motions the radial component of the damping force will contain 



^The discussion above concerns only the lowest-order term in (3, but we shall see below that, to 
all orders, the crucial combination /? can, for circular orbits, be expressed in terms of invariant 
quantities. 
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one power of the small dimensionless quantity R/{Rip). But we learned above, from the 
gauge dependence of the lowest-order damping force, that we can change the definition of 
the radial coordinate so as to set, for instance, the quantity + 1 to zero (for circular 

orbits) . This means that the RHS of Eq. (p.l8|) can be arranged, in the case of quasi-circular 



orbits, to contain three powers of the small parameter R/{R(p). (From Eqs. ( p.ll|) , (|3.12| ) 
we see that for quasi-circular orbits a + (3 oc R^.) We have checked that the reasoning made 
above, using the lowest-order gauge dependence, can be formally extended to all higher PN 
orders. 

The conclusion is that there should exist a special coordinate gauge where, for quasi- 
circular motions, an excellent approximation to the damping force is obtained by replacing 
the radial component simply by zero: 

jr^-^ = Q = jr^--^ _ (3.19) 

To test, a posteriori, the robustness of the approximation (|3.19|) , we shall also consider 
another special gauge: namely that where (a//3)circ = 0. [As we said above, this can be 
achieved at lowest order by a suitable choice of a, and this can be extended to higher PN 
orders by suitable choices of higher gauge parameters.] Finally, this means that there exists 
another coordinate gauge where, to an excellent approximation, the radial damping force is 
given as 

The results in the two gauges are compared and discussed at the end of Sec. 0. 

What is important for the following is that in both gauges (|3.19|) or ( p.20|) , the knowledge 



of the full damping force can be deduced from the sole knowledge of JF 



B. Non-perturbative estimate of the angular momentum reaction force along 

quasi-circular orbits 

The analysis of the previous subsection has shown that the crucial equation in which one 
should accurately incorporate radiation reaction effects is 

dP 

^ = J-|^^-(i?,(^,P«,P^). (3.21) 
12 



As P^p is just the total angular momentum of the binary system, Eq. (|3.21|) expresses the 
rate of loss of angular momentum under gravitational radiation reaction. As usual we shall 
estimate the RHS J^^p = J-'^^^ = T^^^ (remember that JF^*^™ and T^^^ differ only in the 
arguments in which they are expressed) by assuming that there is a balance between the 
mechanical angular momentum lost by the system, and the flux of angular momentum at 
infinity in the form of gravitational waves. In the case of interest here of quasi-circular 
orbits we expect that, to a good approximation, JF^ will not depend explicitly on and will, 
therefore, be expressible in terms of the orbit- averaged finx of angular momentum. Moreover, 
in the case of quasi-circular orbits there is a simple relation between angular-momentum-loss 
and energy- loss. Indeed, the rate of energy-loss along any orbit, in polar coordinates, is given 
by 

f^f^fi^« + ^^„ (3.22) 
and in particular along quasi-circular orbit we have (remembering Eq. ( ^.181 )) 



(^) ^0^f'=+(9(i?^). (3.23) 

\ / quasi— circ. 

Finally, if we know some good estimate of the (averaged) energy-loss along circular orbits, 

say 

we can obtain a good estimate of the needed y^-reactive force 

^circ. ^_£circX^_ (3.25) 

The problem of giving a non-perturbative, re-summed estimate of the energy-loss-rate (or 
"flux function") along circular orbits, say $circ; has been recently tackled by Damour, Iyer 
and Sathyaprakash f^. By combining several of the non-perturbative techniques recalled 
above (to work with an invariant function F{v), to use some global information about F{v) 
in the complex f-plane, to use Fade approximants) Ref. came up with the following 
expression for $circ, considered as a function of the gauge- invariant observable 

v^ = {GMuj/c^y^^; uj = (p, (3.26) 

namely. 
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^ TP r 32 2 10 fDisivu;]^) 



1 - V^/v^o\c{t^) 



(3.27) 



Here, and in the following, we set c = 1 to simplify formulas. The function /dis('^w; ^) 
entering Eq. ( p. 27] ) is the "factored flux function" of , scaled to the Newtonian (quadrupole) 
flux (hence the caret on /dis)- It was shown in [^] that the sequence of near-diagonal Pade 
approximants of /dis(^^) exhibits a very good convergence (at least in the i/ = limit where 
high-order PN expansions are known [^]) toward the exact result (numerically known when 
z/ = [^). On this basis, it was argued in [^] that, in the comparable-mass case, z/ 7^ 0, our 
"best estimate" of / is obtained by Padeeing the currently most complete post-Newtonian 
results, namely the 2.5PN ones p6|. This yields a result of the form 

1 



1 + 



i+- 



Cl V 
ci9 ■)} 



(3.28) 



1+ 



^4 u 



where the dimensionless coefficients q depend only on v. The c^'s are some explicit functions 
of the coefficients fk of the straightforward Taylor expansion of f{v). In turn, the /^'s, being 
defined by the identity (where T means "Taylor expansion" ) 



T [f{v)] ^ T 



1 - 



^pole 



Fiv) 



l + hv + f2V^ + 



are given by 



(3.29) 



fk = Fk- Fk-i/v- 



pole 



(3.30) 



in terms of the Taylor coefficients of the usual (Newton-normalized) flux function 

5G 



T [Fiv)] ^ T 



32 



Fiv) 



1 + F2V' + FsV" + ■ 



[Note that Fi = 0, but that fi = — 1/fpoic 7^ 0.] More explicitly we have 

1247 35 ^ ^ 

F2 = , F3 = 47r , 

336 12 



(3.3i; 



(3.32) 



44711 9271 



65 



9072 



z/ H u , 

504 18 ' 



8191 535 



672 



(3.33) 



and 
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f n h /l /s — /I ro OA\ 

ci = -/i , C2 = /i - — , C3 = -rjj^ TT , (3.34) 

Ji Ji Ui - h) 

/l(/| + /| + /lV4-/2(2/l/3 + /4)) 

^ (/f-/2)(/l/3-/|) ^ J ^^-^^^ 

_ (/i ~ /2) (—/I + 2/2 /a /4 — fi fi — /I /s + /i /a /s) /o oc\ 

(/l/3-/|)(/| + /| + /f/4-/2(2/l/3 + /4)) ■ ^ ^ ^ 

As is clear from these expressions, they depend on the definition used for the quantity fpoio(^) 
which represents a i/-dependent estimate of the location of the "pole" in $circ5 which coincides 
(see the discussion in ||^) with the location of the "light-ring" or last unstable circular orbit 
( ^UgM-ring = 3G'M in the u limit). Actually, as we shall use the Pade representation 
only above and around the LSO (-Rlso'' — 6GM when = 0) the precise choice of Vpoie{i^) is 
probably not crucial (as long as it stays near its known u = limit: fpoic(^ = 0) = 1/v^)- 
In this work, we shall follow Ref. P] and use the pole location they obtained from Padeing 
their "new" energy function e(x), namely 




<^ = ^\l^^.- (3-37) 
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Then, combining Eqs. ( p.25| ), ( p.26| ) and (|3.27|) we define our best estimate of the ip- 
component of the radiation reactive force along quasi- circular orbits as: 

™. _ GM 32 7 his{vu.;jy) . . 



To ease the notation we shall work in the following with reduced quantities, that is: 

_ R — — — ^ — ■ (Q on\ 

1 j-T-improved tt 



^-rM^ , /^eff^^. (3.40) 

G M fj, ji 

Finally, the dynamics, including radiation reaction, in re-scaled coordinates, is explicitly 

described by the following system of equations (in the "canonical" case where Eq. (|3.19| ) 
holds) 

dr dH 

3J = T^O^,Pr,P^ , 3.41 
dt OPr 

dip ^ dH 

-^=^ = -^{r,Pr,P^), (3.42) 
dt op^ 



15 



di 



with 



H 



l + 2v 



'A(r) 1 + 



Bif 



(3.44) 



(3.45) 



T 



32 7 fTiis{.v^]v) 



vv 



(3.46) 



and where in Eq. (|3.45| ) we use the scaled versions of our current best estimate of the effective 
metric coefficients A{r), B{r), see [0 and Eqs. ( |2.14| )-( p.l6| ) above, that is 



A{r) 



2 2iy 



^ ' A(r) 



(3.47) 



Note that the argument entering Eq. ( p.46[) , is simply defined as v^^ = cD^/^, where 
tD = ijj{GM) is the function of r, p^ and p^p defined by Eq. ( ^.42| ), i.e. (2i{r,pr,Pip) = 
dH{r, Pr,P^)/dp^. 



IV. TRANSITION BETWEEN INSPIRAL AND PLUNGE 



The first-order evolution system (|3.41| )- (p.44| ) defines our proposed best estimate for 
completing the usually considered "adiabatic" inspiral evolution into a system which exhibits 
a smooth transition between inspiral and plunge. The rest of this paper will be devoted to 
extracting some of the important information contained in this new evolution system. Before 
coming to grips with such detailed information, it is useful to have a first visual impression 
of the physics contained in our system ( 3.41| )~ (|3.44 ). To do this we plot on the left panel 
of Fig. D the result of a full numerical evolution of Eqs. (|3.41|) - (|3.44|) in the equal-mass case 
{u = 1/4). We started the evolution at r = 15, <y9 = and used as initial values for 
and Pr the ones provided by the adiabatic approximation (see Eqs. ( [4.6|) and (|4.13 ) below). 
The dashed circle in this plot indicates the radial coordinate location of the LSO defined by 
the conservative part of the dynamics, i.e. by the Hamiltonian H{r,pr,Pip). More precisely 
this "r-LSO" is simply defined (for any z/) by r = rLso('^), where rLso('^) is the solution of 
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FIG. 1. On the left panel we show the inspiraling circular (relative) orhit for v = 1/4. The 
location of the r-LSO, defined by the conservative part of the dynamics, is also indicated. On 
the right panel we compare the two kinetic contributions that enter the Hamiltonian: the "radial" 
and the "azimuthal" one. The figure shows that the assumption we made of quasi-circularity, i.e. 
p'^/B{r) <C pt/r"^, is well satisfied throughout the transition from the adiabatic phase to the plunge. 



Eq. (|2.17|) . In particular, tlso {\) = 5.718, as recalled in Eq. (|2.18|) . Note that, in presence 
of radiation reaction effects, there is an arbitrariness in what one would like to mean by 
saying: "the system is crossing the LSO". Indeed, we could define the "LSO-crossing" in 
several inequivalent ways, notably: (i) r-LSO: the time when r = ri^soi'^)] (ii) J-LSO: the 
time when = j = jhsoi^^)', (hi) cu-LSO: the time when dip/di = Q = cDlso('^)- [The 
"LSO" functions of u being defined by solving Eq. ( p.l7|) ; see Eq. (|2.18| ).] This arbitrariness 
is not a problem. Our new evolution system ( 3.41| )-( |3.44 ) describes a smooth transition 
"through" the formally defined "old" LSO, and does not care about old definitions. In 
other words, when u is finite, and especially when ~ 1/4 (which, one should remember, 
is expected to be an accumulation point of observed values of z/; see footnote ^ above) the 
smooth transition process blurs the notion of LSO. It is only for z/ <^ 1 (see below) that one 
recovers a sharp transition near the if-defined LSO. On the right panel of Fig. |I| we compare 
the two kinetic contributions to the Hamiltonian ( |3.45| ): the "azimuthal" contribution p^/r^, 
and the "radial" contribution pi/ B{r). One sees on this Figure that our basic assumption 
of quasi-circularity (which, at the level of H, means pf./B{r) <^ pt/^"^) is well satisfied 
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FIG. 2. Variation with v of the iv-LSO values of the real reduced non-relativistic energy 
-^reai//^ ~ ("^rcai — M) / (on the left), and of the real angular momentum j = P^/ {ixGM) (on 
the right), computed integrating the full dynamics, i.e. with radiation reaction effects included. 



throughout the transition. In fact, even down to r ~ 3.79 one has pf./B{r) < O.lp^/r^. 
The radial kinetic energy would become equal to the azimuthal one only below r = 3. We 
shall, anyway, not use, in the following, our system below the (usual) "light-ring" r ~ 3 
(where pl/B{r) ~ 0.30 p^/r^). 

We exhibit more quantitative results on the transition between the inspiral and the 
plunge in Figs. || and |^. These figures plot the values of several physical quantities (energy, 
angular momentum, radial velocity and radial coordinate) computed at the cj-LSO (i.e. when 
uj = ti;Lso('^)) after integration of the system ( |3.41D -( ^.44D . The energy which is plotted is 
the reduced non-relativistic real energy, i.e. (£^reai — M)//i. [In the test-mass limit, this 
reduced energy equals a/8/9 — 1 = —0.057191.] 

Having obtained, through Figs. ^ and a first impression of the physics of the inspiral 
plunge transition, we shall now study in more detail this transition, notably by comparing 
it with various analytical approximations. The first approximation we shall consider is the 
current standard one used for dealing with the inspiral phase: the adiabatic approximation. 
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FIG. 3. LO-LSO values of the radial velocity (on the left) and of the radial position (on the right) 
versus u, derived integrating the full dynamical evolution. 

A. Comparison with the adiabatic approximation 



Let us compare the exact numerical evolution with the usual adiabatic approximation to 
inspiral motion. This approximation is defined by saying that the (effective) body follows an 
adiabatic sequence of exact circular orbits whose energy is slowly drained out by gravitational 
radiation. It is obtained from Eqs. ( |3.41| ), ( |3.44| ), by neglecting p^, i.e. by setting pr = 0. 
Noticing that dH/dpr = 2pr dH / dpi oc Pr we get that dr/di vanishes linearly with pr. The 
first equation ( |3.41|) is then formally satisfied with pr = = r. Imposing now pr = in 
Eqs. (|3.42| ) and (p.43|) we obtain two further equations: 



dr 



ir,p^) = 0. 



dp. 



[r, p, 



(4.1) 
(4.2) 



where we define 



Ho{r,p^) = H(r,pr = 0,p^ 



\ 



2u 



'A(r) 



1 + ^1-1 



(4.3) 
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Eq. (|4. 1| ) provides a link between r and = j in the adiabatic limit. From the structure 
( ^■45| ) of H, it is easily seen that Eq. ( |4.1| ) is equivalent to looking for the minimum, say (for 
convenience) in the variable m = 1/r, of the "radial potential" 

Wj{u) = A{u) [1 + 3^ u^] . (4.4) 

Solving du Wj{u) = gives a parametric representation of in terms of u: 

"2 / \ (^) / \ 

JadiabN = -^^2^(^)y' (4.5) 

where the prime denote d/du. In the case where the function A is given by Eq. ( |3.47| ), i.e. 
A{u) = 1 — 2u + 2uu^, Eq. ( [4.5|) yields, in term of the orignal (reduced) radial variable 

r = 1/u 

Note that there exist real circular orbits (though possibly unstable ones) as long as Jadiab(^) > 
0, i.e. as long as — 3r^ + 5z/ > 0. In fact the positive, real solution in r of 

[r=^ - Sr^ + 5z/]iight-rmg = (4.7) 

defines the light-ring or last unstable circular orbit (with j^(riight-ring) = +oo). We find 
'^light-ring — 2.84563 in the case u = 1/4. Eq. ( |4.2D then gives the parametric representation 
oiuj = uj (GM) throughout the adiabatic phase for circular orbits: 



1 Vl - 3z//r2 ^ ^ 

Wadiab.(^) = -JH I ' (4-8) 

^l + 2Kv^-l) 

where z{r) denotes the following quantity 

^(^) = H^s{r,Pr =0,p^= Jadiab.) = ^ ^ 2 i " (^'Q) 

r'^ — 3 + 5z/ 

Note that the effective one-body description seems to become somewhat unsatisfactory at 
the light-ring (at least for exactly circular orbits). Indeed, we see from Eqs. ( |4.8| ) and ([4.9| ) 
that the blow up of z{r), i.e. of the effective energy, at the light- ring, Eq. ( [4.7| ), implies 
that the real orbital frequency of circular orbits, uJciTc.{^)y Eq. (^^ ), tends to zero at the 
light-ring. This is probably an unphysical behaviour [from the test-mass limit, one expects 
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the orbital frequency to have a non-zero hmit at the hght-ring; see, e.g., Ref. 0] where 
Pade approximants are used to compute a finite value of ct'iight-ring('^)]- The other factors 
in Eq. ( |4.8| ) imply, as expected, a regular increase of cD(r) as r decreases below the LSO. 
Pending the construction of an improved version of the effective one-body approach which 
would be better behaved, we have decided, when dealing with the evolution of the system 
( p.41| )- (|3.44| ), to stop the simulation at the light-ring. [In our simulations of plunging orbits 
the effective energy stays bounded, but the orbital frequency uj{t) levels off very close to the 
light-ring.] 

Finally Eq. (|3.44| ) becomes in the adiabatic limit 

(4.10) 

Then using dj/di = (dj/dr) (dr/dt) and dip = ujdi we can solve the motion in the adiabatic 
limit by quadratures: 




rfW=( %t) ^-A— . (4.11) 

I J^^(u;adiab.(r)) 



, / 0(7adiabA '^adiab.('^) , I a ^n\ 

(iV'adiab. = I —r- ^— — dr . (4. 12) 

V "'^ / J^<^(a;adiab.(r)) 

The radial velocity Vr = dr/di, as a function of the parameter r, in the adiabatic limit, is 
given by: 

^adiab. ^ .^^(^adiab (r)) 
ajadiab./ctr 

Note that v^'^^^^- formally tends to — oo when r tlso (indeed, jadiab.(^) reaches, by defini- 
tion, a minimum at r = tlso)- This shows that the adiabatic approximation is meaningful 
only during the inspiral phase (i.e. "above" the LSO). In Figs. ^, ^we compare, for u = 1/4, 
the number of gravitational cycles, defined by Mqw = V^gw/(2 vr) = (p/n, the orbital angular 
frequency u (or, equivalently, the gravitational wave frequency, /qw = '^Gw/(27r) = uj/tt), 
and the radial velocity, computed with the exact equations of motion and in the adiabatic 
limit, as well as the gravitational waveform. These figures show that, in the equal-mass case 
u = 1/4, the adiabatic approximation starts to significantly deviate from the exact evolution 
quite before one reaches the LSO. Fig. ^jis normalized so that Aq^^*^ and Af^'^^ coincide for 
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FIG. 4. We compare the number of gravitational cycles (on the left) and the radial velocity 
(on the right), computed with the exact evolution and within the adiahatic approximation, versus 



R/GM. 



large values of R/GM, and that A/"^^^^ be zero at the r-LSO. For instance, we find that the 
number of GW cycles given by the adiabatic approximation differs from the exact number 
already by 0.1 when r ^ 8.8, and that M^^^ (ri^so) = 0.9013 . The left panel of Fig. |^ con- 
trasts uj/u-Lso (= /gw//gwlso where /gwlso = 6"^/^/GM7r is the fiducial Schwarzschild 
LSO GW frequency), computed with the exact evolution and within the adiabatic approx- 
imation, as a function of time. Note that, for the horizontal axis we use cDLSo(0)(t — tLso), 
where ^iLso(O) = T^fGw\so = 6"^^^ (provided by the z/ ^ limit of Eq. (|CT| )) and ^lso is 
defined as the time at which the adiabatic solution reaches the r-LSO position. Finally, on 
the right panel of Fig. ^ we compare the last few GW cycles of the exact and the adiabatic 
restricted waveform, i.e. h(t) = cos 0gw(^)) with v = (dip/diy^^ and 0gw = 2^2, in the 
crucial interesting region around the LSO. By adiabatic restricted waveform we mean the 
restricted waveform in which Lp(i) = adiab. (^) is derived by integrating the two equations 
( [4. lip and ( |4.12| ) (which give a parametric representation of tadiab. (^) and i^adiab. (r) in terms 
of the auxiliary parameter r). 

Note in Fig. ^ that the dephasing between the two waveforms becomes visible somewhat 
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FIG. 5. We contrast the orbital frequency (on the left), divided hy the Schwarzschild value 
^LSo(O) = 6~'^/^, and the restricted waveform (on the right), evaluated with the exact dynamical 
system and within the adiabatic approximation. Note that in both plots the quantities are given 
as a function of the rescaled time variable cDLso(0)(i — ^lso)j where tLSO defined as the time at 
which the adiabatic solution reaches the r-LSO position. 

before the LSO (we shall dwell more on this subject in Section 0). Note also that the time 
when the adiabatic evolution reaches the LSO ("adiabatic LSO") corresponds to a time 
when the exact evolution reaches a frequency u ~ 0.80 c<jlso(0), i.e. a time significantly 
before the tu-LSO. This is why there are more cycles after the adiabatic LSO in Fig. ^ 
(more than two cycles), than there will be after the (exact) cj-LSO (we shall see below that 

jyafterLSO ^ 2ArafterLSO ^ ^ 2048 for U = 1/4). 

B. The r-linearized approximation 

The previous subsection has shown the severe shortcomings of the adiabatic approxima- 
tion. Let us now consider a second analytical approximation which is more accurate than the 
adiabatic one, and which, in particular, allows one to see analytically what happens during 
the transition between the inspiral and the plunge. This approximation is based on a simple 
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linearization with respect to the radial velocity dr/di, which is small during the inspiral, as 
well as the beginning of the plunge. 

As H depends quadratically on pr and 1 we pose 



1 dH 



Pr OPr 



(note that Cr is a positive quantity), where 



uHo{r,j) ifOf,(r,j) 1-6^A" 



(4.14) 



H'^sir,j) = H,s{r,Pr = 0, j) = J A{r) ( 1 + 4 ) • (4-15) 



Then, modulo pf. fractional effects that we neglect, we can write 

dr 



Crir,j)pr. (4.16) 
dt 



Differentiating twice the above equation with respect to time, we obtain 

d'^Pr 1 d^r 



(4.17) 



dp Cr{r,j) dt^ ' 

when neglecting some nonlinear terms oc (dr/di)'^ and (dr/di) (dj/di). On the other hand, 
taking the derivative with respect to time of Eq. (|3.43|) and neglecting fractional corrections 
of Oi^p"^), we end up with 

To get an autonomous system we further approximate j by solving for j in the lowest-order 
approximation to Eq. ( |3.43 ), obtained by neglecting both pr and dpr/dt. In other words. 



j(r) is obtained, as in adiabatic approximation, by solving Eq. ([4.11). Finally, j ~ jadiab.('")5 
as given by Eq. (f4.6|). We define 



u'^ir) = C,.(r, jadiab.('^)) -^^(r, Jadiab.(r)) , 

1 -6r^ + 3i/r3 + 20i/r2 -30i/2 



(4.19) 



Br{r) = Cr(r, JadiabX'")) ■^^(^Jadiab.(r)) J'<^(cDadiab.(?")) , 

2jadiab.(r) [r^ -?,r^ + ^ 

■ ^v[^aA\^hXr)) , (4.20) 



2^^o'(^,Jadiab.) r^r^-Qv) 
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(where the replacements j — > jadiabX'") done after the partial differentiations). It is easily 
seen that the quantity d^HQ/dr dj is negative, so that (JF^ being also negative) the quantity 
Br given by Eq. ( ^.20| ) is positive. 

Combining Eqs. ( [4.17|) and (|4.18| ), we finally derive the following third order differential 
equation in r: 

^+u;,2(r)^^-i?,(r). (4.21) 



We shall often refer to Eq. (|4.21|) as the "linear r-equation" because it was obtained by 
working linearly in the radial velocity f = dr/dt. [Note, however, that this is a third-order 
nonlinear differential equation in r.] It is easily seen that the quantity ^^1{r) defines the 
square of the frequency of the radial oscillations. As seen in Eq. (|4.19|) it is proportional to 
the curvature of the effective radial potential Ho{r,j) determining the radial motion. Above 
the LSO, i.e. when r > rLso('^)) fhe radial potential has a minimum (defining the stable 
circular orbit with angular momentum j) and, therefore, w^ir) is positive. When r = ri^soi^)^ 
the radial potential has an inflection point (see Eq. (|2.17|) ), and, therefore, uj^lr) vanishes. 
When r < rLso('^)) fhe radial potential is concave, and uj'^{r) becomes negative. [See, e.g.. 
Fig. |l] of for a plot of the shape of the radial potential.] 

Within the same approximation used above (i.e., essentially, neglecting terms which are 
fractionally of order p^), we can finally write the angular frequency along our quasi-circular 
orbits as 

dip dHn , , , s 

^^^(r,Jadiab.(r)). (4.22) 

Note that ip is obtained from this equation by a quadrature, once the radial motion r{t) is 
known from the integration of Eq. ( [4.21| ). 



The conceptually interesting feature of the above "r-linearized" approximation is the 
structure of Eq. ( ^4.21| ). The previously considered "adiabatic" approximation corresponds 
to neglecting d^r/di^ in Eq. ( [4.21|) . We now see that this is a good approximation only 
when the characteristic frequency of variation of the radial motion, defined, say, by c^caract. = 
{d^r / dp) / (dr / di) is much smaller than the frequency of radial oscillations cu^ (determined 
by the restoring radial force ensuring the existence of stable circular orbits). As uj"^ tends 
to zero, before changing sign, at the LSO, it is clear that the adiabatic approximation must 
break down somewhat above the LSO. When it breaks down the "inertia term" d^r/dt^ in 
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FIG. 6. Contrast of the number of gravitational cycles (on the left) and the radial velocity (on 

the right), computed with the exact evolution and the linear-r equation, versus R/GM. 



Eq. ( |4.21|) becomes comparable to both the "restoring force" term cu"^ dr / di and the "driving 
force" —Br coming from gravitational radiation damping. 

In Figs. P we compare the number of gravitational cycles and the radial velocity evaluated 
with the exact evolution and the r-linearized equations. We start the evolution at r = 15 
and fix the initial values of dr/di and d^r/di"^ in the "adiabatic approximation" defined by 
neglecting in Eq. ( [4.21|) the "inertia term" d^r/di^ (and then by differentiating again the 
resulting approximate equation). Moreover, we normalize A/gw^"^ be zero at the r-LSO. 
We derive from the exact evolution A/g\v^*('"lso) = —0.04223. The main conclusion drawn 
from Figs. ^ is that the r-linearized approximation is quite good both during the inspiral 
phase and, more importantly, during the transition to the plunge taking place near the LSO. 
This is interesting to know because it shows that the crucial physical effect that is lacking 
in the usually considered adiabatic approximation is the simple "inertia term" d^r/di^ in 
Eq. ( |4.21| ). Note, however, that in order to add this inertia term it is necessary to have in 
hand the Hamiltonian describing at least the slightly non-circular orbits (the normalization 
of Eq. ([4.21|) crucially depends on the knowledge of which depends both on dH/dp^ 
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and on d'^H/dr^). This being said, we do not, however, recommend to use in practice the 
r-hnearized approximation. Indeed, we think that the "exact" system (|3.41|) - (|3.44|) is a 



more accurate description of the evolution of the system because it keeps all the nonlinear 
effects in p"^. Numerically speaking, it is essentially as easy to integrate the "exact" system 
than its f-linearized approximation, so that there would be anyway no practical advantage 
in downgrading the accuracy of the system (|3.41|) -( ^.44|) . However, we shall see next that 
the r-linearized system can be further used to lead to a simple analytical approach to the 
transition to the plunge in the case where u -^l. 

C. The universal p-equation 

Until now we have been considering the general case where the symmetric mass ratio 
u = mim2/ {mi + 7712)^ can be of order of its maximum value z/max = 1/4. As is clear 
from the results above when 4z/ is of order unity the non-adiabatic aspects of radiation 
damping effects become important in an extended region of order A{R/GM) ~ 1, above the 
standard LSO. On the other hand, we expect that when 4z/ <^ 1 the transition between the 
adiabatic inspiral and the plunge will be sharply localized around the standard LSO, defined 
by Eq. (|2.17| ). Indeed, when z/ is a small parameter, the damping force JF^, Eq. ( |3.46|) , being 



proportional to u, can be treated as a perturbatively small quantity in the evolution of the 
system. Consequently, the "driving force" term, —Br, in the r-linearized equation ( |4.21 



contains the small parameter u. It is then clear that all the time derivatives of r (being 
driven by Bj.) will tend to zero with u. If the coefficient u"^ in Eq. ( [4.21| ) never vanishes it is 
easy to see how one would satisfy Eq. (^4.21|) by solving for dr/dt, while considering d^r/dt^ 
as a fractionally small term (to be evaluated by further differentiating dr/di ~ —Br/ujl). In 
that case, one sees that dr/di would be 0{v) (and d^r/di^ = (9(z/^)) as u ^ 0. However, the 
fact that cc'^(r) vanishes when r = rLso(^) shows that the way dr/dt tends to zero with z/, 
near the LSO, is more subtle. Having understood from this reasoning that, when z/ — > 0, the 
interesting transition effects take place very near the LSO, we now turn to a precise analysis 
of this transition. 

A first method for dealing (when u 0) with this transition would be (as just sketched) 



to continue working with the third-order equation (|4.21|) , considered in the immediate neigh- 



bourhood of r = rLso(i^)- However, it is better (in order not to increase the differential order) 
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to go back to the exact system ( |3 .411) -( ^.441) and to approximate it directly when u —>■ and 

Let us see the consequences of the evolution ( |3.41| )-( p.44| ) when r is very near ri^soi^)- 
To do this it is convenient to introduce some notation. Using, as we did in Sec. |IV B| , the 
fact that H depends quadratically on pr and that Pr <^ I, we define: 

LSO 



1 dH 

Pr dpr 



[r,Pr,J) 



(4.23) 



Note that Cl'^'~' is a number, which depends on u |^. In terms of the previous definition 
(|]T|), one has simply C^^^(z/) = C,.(rLso('^), jLso('^))- Explicitly, it reads 



_z/i/o(r,j)^°ff(r,j)(l-6z//r2) 
In the u = limit this simplifies to 



(4.24) 



LSO 



(4.25) 



The point in having introduced the notation ([4.23|) is that Eq. ( |3.41|) reads simply, when one 
is very near the LSO: 

1 dr , ^ 

This allows us to recast Eq. ( |3.43| ) in the form (after neglecting fractional p^ terms on the 
RHS) 



1 (Pr 



dr 



(4.27) 



Here, as above, HQ(r,j) = H(r,pr = 0,pi^ = j). Then we expand the RHS of the above 
equation around the LSO, i.e. we write 



^ As we consider ^ 1, we could further take the limit z/ ^ in all the quantities which have 
a finite limit as = 0. However, in order not to unnecessarily loose accmacy we shall not do so. 
For instance we shall always consider that ri^Q,o[v) is computed for 7^ 0, though we shall see later 
that the direct i^-dependence in rLso(^) (which is 0{v)) is parametrically small compared to the 
width 0(7^2/5 ^) ^Yie radial axis where the transition takes place. 
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(4.28) 



Keeping the first nontrivial terms in the expansion in powers of 5r and 5j (and neglecting 
subleading terms, such as those of order 0{5r5j), 0{{5j)'^) and 0{{5rY)) one obtains 



dr 



{6rY + 



drdj 



(4.29) 



LSO \ / LSO 

Moreover, near the LSO we can write Eq. (|3.44|) as: 

'dHo 



dt dt 



~ J^Jlolso) , with t^LSO 



dj 



(4.30) 
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This yields 



(4.31) 



where tLso is the time at which j{t) = jLso(^)- Let us also define 



ALSO _ ^LSO 



rLSO _ ^LSO 



drdj 



^^(^Lso) • (4.32) 

LSO \ " / LSO 

The quantity 5^^° is the LSO value of the quantity Br{r) introduced in Eq. ( [4.201 ) above. 
The explicit values of these quantities are 



oLSO. ^ 



(^3 _ + 2u) (-210z/ j2 _ eOz/r^ + QOfr^ - 12 f + 6r^ 
2j (r^ - 2r^ + 2u) (r^ - 3r2 + hv) 



J LSO 



In the u = limit they simplify to 



^(^^LSo) 



(4.33) 



LSO 



A,L^°(0) 



1296 ' 



^LSO(^) .^0 



72V3 



1.052 ■ 10~V. 



(4.34) 



Finally, inserting Eq. ( [4.31| ) into Eq. ( |4.29| ), and replacing everything in Eq. ( [4.27| ) yields 
the simple equation 

d^5r 1 



(4.35) 
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This equation can be recast in a universal form by re-scaling the variables 5r and 5t = 
(^-^Lso)- Indeed, posing 

6r = krP, i -ii^so = si = ktT , (4.36) 

with 

it is straightforward to derive the following "universal p-equation" 

ll.ip^^-r. (4.38) 

The explicit values of the scaling coefficients kr and kt are easily derived from our previous 
results. Let us only quote explicitly their u = limit: 

kr{0) = 1.890 , kt{0) = 26.19 u^^^^ . (4.39) 

Note the interesting fractional scalings kr oc z/^/^, kt oc z/~^/^. 

Let us also note the autonomous (time-independent) equation obtained by taking the 
time derivative of Eq. ( [4.38| ) : 



Eq. (|4.40|) could have been directly derived by considering the r-linearized Eq. ( [4.21|) close 



to r = tlso- There is, however, more information in Eq. ( 4.38 ) because its derivation showed 
that r = marks the moment where j{t) = jhsoi^^)- 

The adiabatic approximation is recovered by neglecting in Eq. (^4.38|) the first term on 
the RHS. This gives 

dp\ 1 1 



Padiab. =v^, i-f) =-^= = . (4.41) 

V^^/adiab. Padiab. 

The universal p and p curves and their adiabatic approximations are shown in Fig. ^ We 
have integrated Eq. (|4.38| ) fixing the initial values (for large, negative r) of p and dp/dr in the 
adiabatic limit provided by Eq. ( [4.41| ). We see from Fig. |^ that the adiabatic approximation 



begins to be unacceptably bad when r ~ — 1. From the integration of Eq. ( [4.38|) we get the 
important numerical values: 
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r = 0: p = 0.8339, ^ = -0.8233, (4.42) 

CLT 

p = 0: r = 0.8226, ^ = -1.267. (4.43) 

(IT 

We recall that r = marks the moment where j{t) = jhsoi'^)^ while p = corresponds to 
rit) = tlso- The values given by Eqs. (|4.42|) and (|4.43| ) can then be used to compute corre- 
sponding values of the physical quantities r, dr/dt and j by using the following parametric 
representations derived from our treatment above: 



r(r) = rLso('^) + K p(r) , t(r) = tLSO + ktT, (4.44) 

'dr\ kr dp 

di ) kt dr 



j{r) = jLso(i^) + J'A^LSo) ktT, ( _ ) (r) = ^ ^ . (4.45) 



Correspondingly to these approximate results for r, i, j and dr/di, one can also write an 
approximate result for the angular frequency, namely 

^) M=D(r)=^Wr),,M) 

=.S.soM+f|^) krP(r)+('^) ?l^°k,r. (4,46) 

V / LSO \ / LSO 
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In the approximation where we replace v by zero in all quantities which have a finite limit 
when 0, the above parametric results give the following explicit numerical links [except 
for tLSO which is an arbitrary integration constant] 

r(r) = 6 + 1.890 z/2/^p(r) + 0(z/), t(r) = tLso + 26.19 r , (4.47) 

j(r) = v^- 0.3436 z/^/V + 0(z/), f ^Vr) ^ 0.07216 z/=^/^ ^ , (4.48) 

\dt ) dr 

€j{t) = - 0.03214 p(r) - 0.005062 u^/^ r + 0{u) . (4.49) 
6v6 

Note that these explicit results are less accurate than our previous implicit expressions 
Eqs. ( [4.44| )-( p:.45|) (because of the O^u) error terms entailed by ri^soi^) = ''^lso(O) + O^u), 
etc.). For consistency with the rest of the paper, we have used here (as in Eq. (|4.34|) ) the 
u limit of the value of z/~^ JF^^° defined by the 2.5PN Fade estimate (|3.46| ), namely 
z/~^ J-'^^'~' ~ —0.01312. Note that a more accurate value of this quantity is, according to 
Foisson's numerical results [^] z/^^jF^^° ~ —0.01376, which is ~ 5% larger (in modulus). 
Note the various scalings with u implied (when considering a point in the transition region 
parametrized by some fixed numerical values of p and r) by Eqs. ( |4.44 )-( 4.45 ): notably 
6r = (9(i^^/^), 5j = 0{h'^f^) and pr ^ r = (9(z/^/^). We shall discuss below in more details 
some of these scalings. 

Fig. ^vividly illustrates the fact (mentioned above) that the definition of "LSO-crossing" 
becomes ambiguous in presence of radiation damping. Indeed, for instance, the time where 
= f^Lsoi'^) ("t-LSO"), i.e. the time where p = 0, differs from the time where j = jLso('^) 
("j-LSO"), i.e. the time where r = (see also Eq. ( ^.42| )). An important issue is the 
domain of validity of the universal p-equation, i.e. the range of values of u for which one can 
use Eqs. ( |4.44 )-( ^.45 ) to approximate the transition between inspiral and plunge. We have 
investigated this question numerically by comparing the radial velocity computed with the 
"exact" evolution ( |3.41| )-( |3.44| ), and with the p-equation ( [4.38| ). Let us define the practical 
limit of the domain of validity of the p-approximation by requiring that the fractional error 
in dr/di at the (say) r-LSO be 10%. We find that this limit is reached when u gets as large 
as 

z/„,ax^ 0.05. (4.50) 

Therefore, the explicit expressions above can be used to estimate numerically the physical 
quantities in the transition region only for z/ < i^max- Note that the accuracy of the p-results 
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FIG. 8. Plot of the radial velocity computed both with the exact evolution and with the 
p-equation in the case v = 0.05. The fractional error in dR/dt at the r-LSO is ~ 10%. 



above is, by construction, limited to some small neighbourhood of the LSO. They should 
not be used (even ii u < t'max) to estimate, for instance, the radial velocity at a radius which 
is significantly different from tlso (say at r = 5 or r = 7). To illustrate this we compare in 
Figs. H and ^ the radial velocity computed with the exact evolution, with that deduced from 
the p-equation. We examine three cases: u = z/max = 0.05, z/ = 10~^ and z/ = 10~^. 

Note that, though the accuracy of the approximation defined by the p-equation increases 
as z/ — > 0, its domain of validity actually shrinks as u gets small. Indeed, if we keep p finite 
we see that 5r ~ 1.890 z/^/^ p, tends to zero with u. 

Before discussing the scaling predictions made by the p-approximation, let us comment 
on the various possible definitions of "LSO crossing". We recall that we define: (i) the "r- 
LSO" (by the requirement r(t) = rLso('^))5 (h) the "j-LSO" (j(t) = jLso('^))5 and (iii) the 
"a;-LSO" (cD(t) = a;Lso('^))- [In addition, one can also define an "energy-LSO", and a "naive" 
LSO such that R = 6GM.] We see from our results above that the r-LSO corresponds (in 
the p-approximation) to p = 0, while the j-LSO corresponds to r = 0, and the cj-LSO to 
p -|- 0.1575 z/^/^ r = 0. From these results and the results displayed in Fig. |^ and Eqs. ( |4.42| ) 
and (^4.43|) , we have the following ordering between these LSO's: cj-LSO < r-LSO < j-LSO, 
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FIG. 9. We compare the radial velocity evaluated with the exact dynamical system and with 
the p-equation in the cases v = 0.01 and v = 0.0001. 



where the order symbols refer to the location on the radial axis. We see also that when 
;y2/5 ^ ]^ tu-LSO nearly coincides with the r-LSO. When discussing scaling relations 
it would be essentially equivalent to use any definition of LSO-crossing. For definiteness, 
and for consistency with the rest of this paper where we shall use it, we shall consider the 
cu-LSO (because it is more invariantly defined than the r-LSO). To sufficient approximation 
for determining the leading scaling with u, we shall consider that the a;-LSO corresponds to 
p ~ 0. 

One of the most useful scaling law to consider is that concerning the radial momentum 
at the tu-LSO. Combining Eqs. (|4.26|) and ( [4.45| ) we get: 



dr 



Pr 



^ ( /iLSOn-2/5 / rLSOn3/5 '^P 



(4.51) 



From Eq. (|]4|) the value of dp/dr at the w-LSO (i.e. p ~ 0) is dp/dr ~ -0.8233. Using 
also the numerical values (taken when z/ — 0) of the coefficients entering Eq. (|4.51|) , we get 
the predicted scaling 



{Pr)u-i.so ^ -0.0844 (4z/) 



3/5 



(4.52) 
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FIG. 10. Scaling laws for the radial momentum and the post-LSO number of orbits provided 
by the p-approximation. On the left panel we show the exact numerical results for the radial 
momentum obtained by integrating the full evolution system down to the iv-LSO. Two fits of the 
type pr = —a (Av)^^^ and pr = —a (4z^)* are also indicated. On the right panel the number of orbits 
remaining "after LSO-crossing" is compared with the numerical results computed from the exact 
evolution. We have indicated both the best fit to a formula of the type A(p/2'iT = a (4i^)~^/^, and 
of the type Aip/27r = a(4z^)^. Note that, even if the figure covers the range of values of u up to 
1/4, both fits have been evaluated including values only up to z^max = 0.05. 



In the left panel of Fig. we compare the analytical scaling prediction, (pr)w-Lso oc (4z/)^/^, 
with the numerical results obtained by integrating the full evolution system ( p.41| )- (|3.44| ) 
down to the cj-LSO. We have also computed the best fits to the data using either a formula 
with one free parameter, of the type Pr = —a (4z/)^/^ or with two free parameters, Pr = 
—a (4z/)*. Note that the predicted scaling is a surprisingly good fit to the exact results, even 
for values of u much larger than the domain of validity of the p-equation. In fact, it is 
numerically quite accurate even for z/ = 1/4. [In the one-parameter fit, note that the best-fit 
coefficient a = 0.0750 is 11% smaller than the calculated one, Eq. ( [4. 521 ). This is because the 
best-fit one takes into account the values of Pr for larger values of u than the test-mass-limit 
result iH).] 
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Another useful scaling law concerns the number of orbits remaining "after LSO-crossing" . 
Let us define the number of orbits after LSO-crossing as A(/9/27r where Ay9 is the difference 
in orbital phase between the "light-ring" r = rught-ring (obtained from Eq. (|4.7])) and the 
c<j-LSO, u = ui^soi^)- This quantity cannot be really estimated within the p-approximation, 
because this approximation assumes that 6r <^ 1. However, we can formally say that, 
within the p-approximation, we wish to consider the asymptotic limit where p tends to 
— oo proportionally to z/^^/^ (so that 5r is finite). The question is therefore: what is the 
asymptotic behaviour of the solution p = pij) of Eq. ( [4.38| ) when p — oo? It seems that in 
this limit the "source term" — r on the RHS of Eq. ([4.38|) is relatively negligible. Indeed, let 



us neglect it and solve the approximate equation f;^ + | = 0. This equation describes the 
motion of a particle (p = —dV{p)/dp) with potential energy V{p) = p^/6. This potential 
energy (which represents the effective radial potential near the inflection point corresponding 
to the LSO) is unboundedly negative when p — oo. Writing the conservation of "energy", 
I p^ + V{p) = const, one finds that, as p ^ — oo, the kinetic energy grows without bound 
and approximately satisfy ^ p^ ~ ~^(p) whose solution is 

p = -12 (t^-t)-^ (4.53) 

for some constant t^o- We conclude that, as p — oo, the variable r tends to a finite limit 
Too- [We find Too — 3.9. The corresponding curve is shown in the left panel of Fig. |^.] 
Therefore, from Eq. (|4.53| ), the total time elapsed after the LSO, too — ^lso, scales like z/^^/^. 
Correspondingly, within the p-approximation, the leading approximation to the orbital phase 
(obtained by integrating the zeroth order term in Eq. ( 4.53|) ) reads 



As 

cDlso admits a finite limit as z/ — > 0, we expect from Eq. ( [4.54| ) the scaling law 



^oc(4z/)-^/^ (4.55) 

ZTT 

This prediction is compared in Fig. |10| with the numerical results obtained by integrating 
the full system ( |3.41|) -( ^.44|) . As expected from the necessity to inconsistently consider 
parametrically large values of p oc z/^^/^, this prediction is less accurate than that obtained 
for the radial momentum at the a;-LSO. We have indicated both the best fit to a formula of 
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the type /\ip/2Tx = a (4z/)"^/^, and the best fit to A(p/2'K = a (4z/)^. Note that, both fits have 
been evaluated including values of u only up to z/max = 0.05. Indeed, as discussed above, 
beyond this value the fractional error in the radial velocity at the tlso is ~ 10%. 

Some comments are in order concerning these results. First, we note that although 
A^Lso ~ A(f/27T tends to infinity when u ^ 0, it does so very slowly so that the total 
number of orbits after the LSO is always quite small compared to the number of orbits "just 
before and around the LSO". Let us define the latter number as NI^q^'^ = /orbit. //orbit. = 
i /gw//gw where /orbit. = | few = uj/2tt denotes the orbital frequency, and /orbit, the time 
derivative of the orbital frequency caused by GW damping. In the adiabatic approximation, 
combined with a Newtonian approximation for both the orbital energy and the GW fiux, 
this number reads (see, e.g., 0) 

2 924 

N^sT"^-^- (4.56) 

The ratio A/'£s67^lsT'^ - 0.3446 (4z/)^/^ (derived using the result of the fit, i.e. A^gfg' = 
1.0075 (4z/)^^/^) is therefore parametrically small as z/ — 0. This suggests that, when 4z/ <^ 1, 
the existence of even a formally parametrically large (oc u'^^^) absolute number of cycles 
left after the LSO will have only a fractionally negligible effect on the extraction of a GW 
signal from the noise by means of relativistic filters built on the adiabatic approximation, 
and terminated at the LSO H, If. On the other hand, when 4z/ ~ 1 the ratio A^^fo /iV^so"*^ 
is not very small. In particular, when z/ = 1/4 the number of orbits after the w-LSO is equal 
to A/'lso = 1/4) = 0.6024 (computed from the exact evolution), while iVf^so°'^(z/ = 1/4) = 
2.924. The ratio between the two is iV^so/A^f^so"'^ = 0.2060. As recently emphasized in 
Ref. P], the fact that A^f^g'^'^'^ is not large means that the filtering of such a signal out of the 
noise is a delicate matter which sensitively depends on the modeling of the phase evolution 
near the LSO, and on the modeling of what happens to the signal after LSO crossing. In 
Ref. it was assumed that the signal is abruptly terminated at the LSO. In a later section 
we shall use the tools introduced here to go beyond such an approximation and study the 
part of the waveform which is emitted after LSO crossing. 
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V. INITIAL DATA FOR NUMERICAL RELATIVITY 



One of the main aims of this paper is to use the improved approach to the transition 
from the inspiral to the plunge introduced above to compute initial dynamical data (i.e. 
initial positions and momenta) for binary black holes that have just started their plunge 
motion. Ideally, we wish to give dynamical data for two black holes {Qi, q2,Pi,P2) such 
that the coordinate distance \qi — is: (i) large enough that one can trust the re-summed 
non-perturbative technique allowing one to compute these data; (ii) large enough to allow 
one to hope to complete the present work by constructing the initial gravitational data 
{gij{x), Kij{x)) determined (in principle) by {q^,p^)] and, finally (iii) small enough to leave 
only less than one orbit (at least when ~ 1/4) to evolve by means of a full 3D numerical 
relativity code. We think that the point (i) is satisfied if we use the Pade-type plus 
effective-one-body methods we have combined above, and if we stop the evolution of quasi- 
circular orbits anywhere around the LSO. We shall leave the point (ii), i.e. the important task 
of completing the present work by constructing gravitational data to future work. However, 
in preparation for this task we shall show how one can compute the dynamical data {q^, p„) in 
the convenient ADM coordinates. Indeed, the coordinate conditions introduced by Arnowitt, 



Deser and Misner |^ have the double advantage: (a) to be linked to the 3 + 1 formulation 
which is used in numerical relativity, and (b) to be linked to explicit, high-order post- 
Newtonian calculations Jl^. Concerning the point (iii), the work above shows that if we 
stop the inspiral + plunge evolution at the (invariantly defined) tu-LSO (i.e. when dip/dt = 
uji^soi'^)) there indeed remains (when 4z/ ~ 1) less than one orbit to go before reaching the 
light-ring (see next section for a discussion of the importance of the light-ring). Note that 
there is nothing sacred about giving data precisely at the tu-LSO. Because of the points (i) 
and (ii) above we wish to stay "as high as possible". Because of point (iii) we must, however, 
be just after LSO crossing. As was already discussed, there are several possible definitions of 
"LSO crossing". The cu-LSO is the innermost LSO (see below) and is therefore a convenient 
choice (however, there would be nothing wrong in giving data at a slightly different place; 
in fact we recommend to do it to check the robustness of the numerical spacetimes evolved 
from our data). 

As just recalled we wish to (numerically) compute complete dynamical data at the u- 
LSO, and in ADM coordinates. The evolution system ( |3.41 )-( ^.44 ) given above allows one to 
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compute dynamical data (r, (p,pr,p^) for the relative motion described in (reduced) effective 
coordinates (i.e. the coordinates used in the effective-one-body description). In Ref. ^ we 
have shown how to map the ADM positions and momenta (gADM^^ADM-j qj-^^q ^^j^g effective 
positions and momenta {q,p) by means of a generating function G{q^^^,p). Let us first 
recall, in order to avoid any confusion, the trivial transformations linking Cartesian-like to 
polar-like coordinates, as well as those linking the original to the scaled coordinates. We 
recall that we work in the center of mass frame and that we consider planar motion in the 
equatorial plane 6 = 7i/2: 

Q' = <A-(A. Pi= Pii = -P2i , (5.1) 



R 



n'Pi, P^ = Q^Py-Qyp,, (5.2) 



^ Pi = —, (5.3) 



CM' A* 
P P 

Pr = — =n'pi, Pv = —77^ = Py - Px ■ (5.4) 

H jJL (jlVl 

Here it = Q^/R = q^/r is the radial unit vector {R = \Q\, r = \q\). We have also = 
Rcosip, Q'^ = Rsimp, q^ = rcosip, q^ = rsimp. The relations above hold both in effective 
coordinates (denoted by {q\Pi) without extra labels) and in ADM coordinates (i?ADM;Pf^^)- 
The link between {q\pi) and (o'ADM'Pf^'^) defined by a generating function G(g'ADM'P*) 
and reads 

,.,,,,,,^.^(^, (5.5) 

(5,, 

The generating function G has been derived up to 2PN order in [|^ (see Ref. [jl5| for the 
determination of G at the 3PN level) 

G(g^°^,p) = 4 G'iPN(g™,p) + \ G2PN(g™,p) • (5.7) 



The partial derivatives needed in Eqs. (|0| ), (|5.6| ) read 
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dGipn{q,p) 
dpi 



Pt 
q' 



-2'' +0 + 2) 



V\ 1 



p [q-p)'^: 



(5.8) 
(5.9) 



dG2PN{q,p) 
dq' 



Pi 



-z/(l + 3z/)p^ + -(2-5z/)^^ + -z/| 



+ 3z/) 



+ -(l-7z/+z/" 
4 

P 



+ qi{q- p) 



3 

U I 



+ 3z/) 



19 PJ 



^G'2PN(g,p) 
dpi 



z/ (1 + 3 z/) + - (2 - 5 z/) — + - z/ (8 + 3 z/ 



19 PJ 

g3 



+ -(1 -7z/+z/' 
4 ^ 



+ (q ■ p) 



|(l + 3z.)p^ + ^(2-5z.)i 



(5.10) 



(5.11) 

-dH), 



Given and pi, we use first Eq. ( |5.5|) , and the values of tlie partial derivatives (|5J 
to solve numerically for ^Xdm- Then we use Eq. (|5.6| ) to compute 

The initial data we start with are the results of the numerical integration of the system 
( p.41| )- (|3.44| ), i.e. the values of r, yj, pr and p^p at some time in the evolution (which we 
choose to be the time when uj{t) = ui^soi'^))- Actually, the value of (p is without significance 
and we renormalize it to the convenient value v^new = so that we work with Cartesian-like 
data of the simple form (remember that we work in the x — y plane, q^ = = Pz, and that 
we simplify the writing by denoting = when working in Cartesian-like coordinates) 







Px Pr ; 



Py 



P^ 

r 



(5.12) 



ADM ^ADM 



P 



When solving, as indicated above, Eqs. (|5.5|), ( |5.6| ) to derive g^^^, g^^^ and p: 
we get these quantities in a not optimally oriented coordinate system (i.e. though we started 
with qy = 0, we end up with g^°^ 7^ because there is a rotation between the two coordinate 
systems). As the global orientation is of no physical significance, it is convenient to turn the 



ADM coordinate system by an angle a so that ip 



ADM 



ADM, 
'f'old 



a 



0. In other words, after 



this rotation one has, as in Eq. ( |5.12| ) above 



ADM new 

qx 



ADM 



ADM new 



n „ ADM new 
5 Px 



P 



ADM 



ADM new 
Py 



P\ 



ADM 



r 



ADM ■ 



(5.13) 



The angle of rotation a is determined by 
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„ADM 



p. 



ADM 



n 



ADM 



ADM 
Pip 



a 



0.25 

0.1 

0.01 

0.001 
0.0001 



4.717 
4.853 
4.938 
4.948 
4.949 



-0.07570 
-0.04425 
-0.01163 
-0.002992 
-0.0007592 



0.7021 
0.6997 
0.6996 
0.6998 
0.6999 



3.312 
3.396 
3.455 
3.463 
3.464 



-0.006256 
-0.001524 
-4.088 • 10-5 
-1.054 • 10"^ 
-2.675 • 10'^ 



TABLE I. Initial data in ADM coordinates at u>-LSO for five representative values of u. 



,ADM old 
1. 

while the more invariant quantities r^^^ and p^^^ are given by 



^ADM 



/l'„ADM\2 I |'„ADMN2 „ADM ^ / „ADM ADM , ADM ADM\ 

Y \1xold ) "T \%old ) ^ Pr — ^ADM \^xo\d Pxold "r 'iyold Pyold ) ■ [O.IO) 

Note that (because of the rotational invariance of G) all the angular momenta coincide: 

„ _ „ADM _ „ „ _ „ „ _ „ADM ADM _ ADM ADM _ ADM ADM _ ADM ADM /p- . 
Pip ~ Pip — Hx Py HyPx—Hxo\d Pyold 'iyold Fx old ~ Hx new Py new Hynew Pxnew ■ l^^-J-U; 

This relation is a useful check on the numerical precision of the solution of Eqs. ( ^.5|) , (|5.6|) . 

In Tab. | we give initial data in ADM coordinates at the tu-LSO for five values of the 
parameter v. We give the more invariant quantities corresponding to the "new" ADM 
coordinate system Eq. ( |5.13| ). The quantity p^^^ denotes the "transverse" momentum, i.e. 
simply p^DM ^ pADM^^ADM ^ Py^^. For completeness, we give also the value of the angle 
a, Eq. ( CT )- 



So far all the results we have discussed considered the evolution system (|3.41|) - (|3.44|) as 



the "exact" description of the transition through the LSO. However, as discussed in Sec. Ill 
this system is more like a convenient fiducial system within a class of systems obtained by 
shifting (by 0{v^ /c') terms) the coordinate system. To test the robustness of our predictions 
for physical quantities at the LSO we shall now compare the results of the fiducial system 
( p.41| )- (|3.44| ) with the results obtained by the more general system ( ^^ -( ^^ , with a radial 



force given (in terms of JF^^) by Eq. ( |3.2CI| ). For simplicity, we consider only the (most 
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crucial) equal-mass case, v = 1/4. We find that our fiducial system (with jF^j = 0) yields 
the following numerical values at the tu-LSO (when starting with an orbital phase = at 
r = 15) 

r = 5.639, = -0.07432, f = -0.03563, (5.17) 

pNR 

ip = 82.72 , j = 3.312 , = -0.01640 . (5.18) 

On the other hand, the system including the non-zero radial force ( |3.20|) yields at the cj-LSO 
(still starting with an orbital phase = at r = 15) 

r = 5.638, = -0.07388, f = -0.03542, (5.19) 

pNR 

(^ = 82.77, j = 3.311, ^ = -0.01643. (5.20) 

As we see the differences in the numerical results are quite small. For instance, the fractional 
change in the (crucial) radial momentum is less than 6 x 10~^. We note also that the 
dephasing at the LSO is only 0.05 radians. This analysis indicates that the results based 
on our fiducial system are quite robust, mainly because our basic assumption of "quasi- 
circularity" (i? -C i?0) is well satisfied during the transition to the plunge. 



VI. GRAVITATIONAL WAVE-FORMS FROM INSPIRAL TO RING-DOWN 

In this section, we provide, for data analysis purposes, an estimate of the complete wave- 
form emitted by the coalescence of two black holes (with negligible spins). This estimate 
will be less accurate than our results above because we shall extend the integration of our 
basic system ( |3.41|) -( ^.44|) beyond its range of validity. We think, however, that even a rough 
estimate of the complete waveform (exhibiting the way the inspiral waveform smoothly trans- 
forms itself in a "plunge waveform" and then into a "merger plus ring-down" waveform) is 
a very valuable information for designing and testing effectual gravitational wave templates. 
[See, in particular, the recent work which emphasizes the importance of the details of the 
transition to the plunge for the construction of faithful GW templates for massive binaries.] 

Our (rough) assumptions in this section will be the following: (i) we use the basic evolu- 
tion system (|3.41| )- (|3.44|) to describe the dynamics of the binary system from deep into the 
inspiral phase (say r ^ 15) down to the "light-ring" r = riight-ring('^) — 3; (ii) we estimate 
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the waveform emitted during the inspiral and the plunge by means of the usual "restricted 
waveform" approximation 

i < tend ■ /iinspiral(il = C v'^{i) COs(0Gw(^)) , = [ ] ' ^GW = , (6.1) 



, dt , 

and (iii) we estimate the waveform emitted during the coalescence and ring-down by match- 
ing, at a time t = tend where the light-ring is crossed, the inspiral + plunge waveform ( |6.1| ) 
to the least-damped quasi-normal mode of a Kerr black hole with mass and spin equal to 
the total energy and angular momentum of the plunging binary (at t = tend): 

i > tend : /^merger (t) = A e-^*-*-^)/- cOs(o^qnm (t - tend) + B) . (6.2) 



For convenience, we shall normalize the waveform by taking C = 1 in Eq. ( |6.1|) . The 
amplitude A and the phase B of the merger waveform ( p.2| ) are then determined by requiring 
the continuity of h(i) and dh/di at the matching point t = tend. 

Before giving technical details let us comment on our assumptions (i)-(iii). First, we 
recall that Fig. |l| had shown that the quasi-circularity condition pi/ B{r) <^ p'^/r'^ (which is 
the basic condition determining the validity of our evolution system) was satisfied with good 
accuracy during the inspiral and the beginning of the plunge, and was still satisfied, though 
with less accuracy {p1/B < O.Sp^/r^ in the worst case v = 1/4) down to the light-ring 
r ~ 3. In other words, our work is showing that the so called "plunge" following the inspiral 
phase is better thought of as being still a quasi-circular inspiral motion, even down to the 
light-ring. We therefore expect that the usual restricted waveform ( |6.1D (valid for circular 
motion) will be an acceptable description of the GW emission during the plunge. Note 
that we consider that the description of the amplitude of the gravitational wave in terms of 

= ip'^^^, being simpler and more invariant, has a better chance of being correct than a 
description in terms of some other Newtonian-like approximation to the "squared velocity" 
such as (r0)^ or 1/r. Some evidence for this faith is given by the fact that the GW fiux 
is surprisingly well approximated (within 10% down to the LSO) by the usual "quadrupole 
formula" if the velocity used to define the quadrupole formula is the invariant = 0^/^ 
(see, e.g.. Fig. 3 of [|]). 

Concerning the choice of the light-ring for shifting the description between a (quasi- 
circular) binary motion and a deformed Kerr black hole, our motivation is twofold. First, 



in the test-mass limit, ^ 1, it has been realized long ago, in the first work p8| which 
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found the existence of a merger signal of the type (|6.2| ) following a plunge event, that the 
basic physical reason underlying the presence of a "universal" merger signal was that when a 
test particle falls below R ~ 3GM, the GW it generates is strongly filtered by the potential 
barrier, centered around R ^ 3GM, describing the radial propagation of gravitational waves. 
It was then realized [29] that the peaking of the potential barrier around R ~ SGM is itself 
linked to the presence of an unstable "light storage ring" (i.e. an unstable circular orbit for 
massless particles) precisely at i? = -Rught-right = 3GM. A second argument (applying now 
in the equal-mass case, z/ = 1/4) indicating that riight-right(l/4) — 2.84563 is an acceptable 
divide between the two-body and the perturbed-black-hole descriptions comes from the works 
on the, so called, "close limit approximation" [^. Indeed, recent work (see the review 
|?T[) suggests a matching between the two-body and the perturbed-black-hole descriptions 



when the distance modulus jiQ ~ 2. Using the formulas of Ref. one finds that ~ 
2 corresponds to a coordinate distance in isotropic coordinates of r^^° ~ \/2xGM. This 
corresponds to a Schwarzschild-like radial distance ~ r (1 + GM j'lrf' — 2.59GM which 
is not very far from -Riight-ring(l/4) ~ 2.84GM. 

In keeping with our prescription of setting the divide between a binary-black-hole de- 
scription and a perturbed-single-black-hole one, at the time tend, when r ~ Tiight-ring('^), it 
is natural to assume that the final hole formed by the merger is a Kerr hole with mass Mbh 
and angular momentum JTbh given by: 



^ =ifend = -Vl + 2^W -1), ^-d^;:^' (6.3) 

while the dimensionless rotation parameter a is: 

''-=GMi;- 1 + 2.(5- -!)• * ' 

As the system reaches the stationary Kerr state, the non-linear dynamics of the merger 
becomes more and more describable in terms of oscillations of the black hole quasi-normal 
modes [Q. During this phase, often called the ring-down phase, the gravitational signal will 



be a superposition of exponentially damped sinusoids. The gravitational waveform will be 
dominated by the / = 2, m = 2 quasi-normal mode, which is the most slowly damped mode. 

As a rough approximation we assume that the full merger -|- ring-down signal (starting 
when the light-ring is reached) can be represented in terms of this least damped quasi-normal 
mode. If ujqnm denotes the circular frequency of this mode, and r its damping time, this leads 
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FIG. 11. We compare the inspiral + plunge waveform, terminated at the light-ring, to the 
adiabatic waveform, terminated at the adiabatic LSO. 



to the simple description (|6^). The quantities {ugnm^T) are functions of (Mbh^o^bh) which 
have been investigated numerically 0], [Q. Using analytic fits the following expressions 



for the frequency and the decay time of the quasi-normal modes were obtained |^ 

Mbh c^qnm = [l " 0.63 (1 - a)=^/^°] //(a) , (6.5) 

TUJqnm = 4 [1 - a]"^^^° fqCa) , (6.6) 



where //(a) and fQ{a) are correction factors provided by Tab. 2 of I^Sf. Note that // = 0.9587 
and /q = 1.0501 for a = 10"^. 

We have numerically studied only the equal- mass case p = 1/4. We have chosen the 
matching point tend such that r (tend) = night-ring(l/4) = 2.84563. With this value of tend we 
obtain the following values for the characteristics of the formed black hole: 

aBH = 0.7952, Ebh = 0.9761 M, (6.7) 
Mu^nm = 0.5976 , M/r = 0.07795 . (6.8) 

Note the numerical value of the quasi-normal mode frequency 



/,..,„.^^ 1885 ^-^lHz. (6.9) 
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FIG. 12. Plot of the complete waveform: inspiral and plunge followed by merger and ring-down. 
The locations of several possible definitions of LSO crossing are also indicated. 



Our results for the waveform are shown in Figs. |11| and |I2|. In Fig. ^ we compare the 
inspiral + plunge waveform ( |6.1| ) (terminated at the light-ring) to the usually considered 
adiabatic waveform (terminated at the "adiabatic LSO"). As already discussed in Section 
1^, by "adiabatic waveform" we mean a restricted waveform (|6.1| ) (with C = 1) in which 
ip(t) = ipa.dmh.(t) is defined by integrating the two equations ( [4.11| ) and ( [4.12|) . This Figure 
shows that there is a significant dephasing of the adiabatic waveform with respect to the 
(more) exact one already before the LSO. Moreover, the real inspiral signal continues to 
increase and oscillate for ^ 2.35 cycles after the adiabatic LSO. 

In Fig. |l^ we plot our estimate of the complete waveform: inspiral and plunge (solid line) 
followed by merger and ring-down (dashed line). We also indicated the locations of several 
possible definitions of LSO crossing (see Section |V| above). In addition to the definitions 
mentioned above we also included a "naive LSO" (defined simply by r 



= 6 as in the 



Schwarzschild geometry) and an energy-LSO (such that £^rcai(^) = ^reai ('^))- 
The corresponding numerical values of the reduced radial coordinate r are: 



rj.Lso = 6.631 , rllZ = 6.000 , r^_Lso = 6.534 , 
ri._Lso = 5.718 , r^-Lso = 5.639 . 



(6.10) 
(6.11) 



As mentioned above, the fact that the various definitions of the LSO differ significantly is 
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due to the fact that when z/ = 1/4 the GW damping effects are rather large and blur the 
transition to the plunge. Note that the number of GW cycles left after the (exact) tu-LSO 
(and until the light-ring) is A^^^^'' = 2A^„f{;'?[ = 1.2048 (for z/ = 1/4). As said above, this 
is smaller than the (physically less relevant) number of cycles left after the adiabatic LSO 
(where u ~ CSOculso); which is ~ 2.35. 

Even if our estimate of the waveform is admittedly rough, we think that it can play an 
important role for defining better filters for the search of signals in LIGO and VIRGO. In 
particular, two features of this waveform are striking: (i) the 'plunge' part of the waveform 
looks like a continuation of the inspiral part (this is because the orbital motion remains 
in fact quasi-circular), and (ii) the adiabatic waveform gets significantly out of phase with 
the exact waveform before crossing the LSO. We shall come back in future work to the 
consequences of these results for data analysis, and see how they can be used to improve 
upon the state-of-the-art filters constructed in Ref. 0, [§]. 



VII. DISCUSSION 



In this paper we have extended a methodology introduced in previous papers 0,0) and 
applied it to the study of the transition from inspiral to plunge in coalescing binary black 
holes with comparable masses, moving on quasi-circular orbits. Our philosophy is that it is 
possible to use suitably re-summed versions of post-Newtonian results to write an exphcit 
(analytical) system of ordinary differential equations describing the transition to the plunge. 
Our explicit proposal is the evolution system ( p.41|) -( ^.44|) obtained by combining the results 
of for the re-summation of the gravitational wave damping, and the results of [0] for the 
re-summation of the conservative part of the dynamics of comparable-mass binaries. The 
basic reason why we think the simple evolution system (|3.41|) - (|3.44|) can accurately describe 
the transition to the plunge is that we have consistently checked that most of the "plunge" 
motion (at least down to -R ~ SGM) is in fact very much like a quasi-circular inspiral motion 
(with i?2 ^ {R^f). 

In general one needs to numerically integrate the basic evolution system ( p.41| )-( p.44| ) to 
get physical results of direct interest. However, we have shown that one can understand the 
various physical elements entering this system by comparing it to several simple approxima- 
tions: the adiabatic approximation, the r-linearized one, and the universal p-approximation 
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(valid when z/<0.05). In particular, the latter approximation allowed us to derive some 
scaling laws: one scaling law (which is very well satisfied, even up to the maximum value 
V = 1/4) states that the radial momentum at the Last Stable Orbit (LSO) scales like z/^/^, 
while another scaling law (accurately satisfied only for z/ ^ 1) states that the number of 
cycles left after the LSO scales like z/~^/^. 

The two most important consequences of the present approach are: (i) a way to com- 
pute initial dynamical data (qi, Pi? P2) ^ comparable-mass binary black hole system, 
represented in ADM coordinates, such that only a fraction of an orbit needs to be further 
evolved by numerical relativity techniques, and (ii) an estimate of the complete waveform 
emitted by a binary black hole coalescence, smoothly combining an inspiral signal, a plunge 
signal, a merger signal and a ring-down. 

However, much work remains to be done to firm up and complete our approach. We 
checked the robustness of our approach by considering an as-well-justified, slightly different 
evolution system. But stronger checks are called for. In particular it would be quite im- 
portant to extend the present work (which used as input the 2.5PN-accurate damping and 
2PN-accurate dynamics) to higher PN levels, when they become fully available. We note 
in this respect the recent work which extended the effective-one-body approach to the 
3PN level. [Note in passing that quasi-static tidal interactions between black holes enter 
only at the 5PN level [jl2|-] It is quite important to complete our determination of initial dy- 
namical daXa, {qaiPa) by explicitly constructing the initial gravitational data, {gij{x), Kij{x)) 
corresponding to {qa,Pa) (and containing no free incoming radiation). When this becomes 
available it will be possible to further check our method (by numerically evolving space- 
times starting at various stages of the plunge) and to provide more accurate estimates of the 
merger waveform. Though our "light-ring-matching" approach to estimating the complete 
waveform is admittedly rough, we think it can play a useful role for data analysis: it can be 
used to test the accuracy of present templates (based on the adiabatic approximation) and 
allow one to construct more accurate, or at least, more robust, templates. We will come back 
to this issue in future work. Finally, let us note that it would be, in principle, important to 
be able to extend our approach to black holes having significant intrinsic spins. We, however, 
anticipate that this is a highly non-trivial task. 
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